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We have used Monte Carlo simulation techniques to obtain 
the magnetic phase diagram of the double exchange Hamil- 
tonian. We have found that the Berry phase of the hopping 
amplitude has a negligible effect in the value of the magnetic 
critical temperature. To avoid finite size problems in our sim- 
ulations we have also developed an approximated expression 
for the double exchange energy. This allows us to obtain the 
critical temperature for the ferromagnetic to paramagnetic 
transition more accurately. In our calculations we do not 
observe any strange behavior in the kinetic energy, chemical 
potential or electron density of states near the magnetic crit- 
ical temperature. Therefore, we conclude that other effects, 
not included in the double exchange Hamiltonian, are needed 
in order to understand the metal-insulator transition which 
occurs in the manganites. 

PACS number 71.10.-w, 75.10.-b 



I. INTRODUCTION 

The recent discovery of colossal magnetoresistance 
(CMR) m in mixed-valence compounds of the form 
Lal± x Af?Mnlt x Mni + Ol- , (where A can be Ca, Sr or 
Ba) has revived the interest in these materials with per- 
ovskite structure. The interest is focused on the phase 
diagram and the magneto-transport properties. Both 
features have been determined experimentally For 
0.1 < x < 0.5 and low temperatures the system is metal- 
lic and presents ferromagnetic order. As the tempera- 
ture increases the system becomes insulator and param- 
agnetic. The ferro-paramagnetic transition occurs at an 
x dependent critical temperature T c {x) ~ 300X. The 
colossal magnetoresistance effects are observed for tem- 
peratures close to T c (x). For x < 0.1 and low tempera- 
tures the system is a 'layer' antiferromagnetic with fer- 
romagnetic coupling inside planes. 

In these compounds the electronically active orbitals 
are the Mn d orbitals. The number of d electrons per Mn 
is x + 3, namely, four electrons per Mn 3+ and three per 
Mn A+ . The cubic symmetry and the strong Hund's rule 
coupling make that three electrons get trapped in the ti g 
states and therefore these electrons become electrically 
inert, forming a core spin S of magnitude 3/2. The rest 
of the electrons go to the e g orbitals and they get strongly 
coupled to S by the Hund's rule coupling. Consequently, 
all the spins on each Mn prefer to be parallel. For small 
values of x, the perovskites show a long range Jahn- Teller 
order which selects a preferred d orbital and therefore it is 
possible to assume that the electrons move only through 



one d orbital. 

The double exchange (DE) mechanism developed by 
Zener , explains the existence of ferromagnetism and 
metallic behavior at low temperatures. In this model, the 
electrons get mobility between the manganese ions using 
oxygen, which is magnetically inert, as an intermediate. 
As a consequence, the tunneling takes place between two 
configurations in which the Mn ions of different charge 
(Mn 3+ and Mn A+ ) interchange their valence states. This 
conduction process is proportional to the electron trans- 
fer integral t if the core spins of Mn 3+ and Mn A+ are 
aligned. Otherwise, the transfer integral is inversely pro- 
portional to the Hund's rule coupling energy, which for 
the Mn ions is much larger than t. In the DE model the 
ferromagnetism is then induced via this electron conduc- 
tion process. 

In the limit of large Hund's coupling the spin of the 
electron in the active d orbital is parallel to the local core 
spin Si, which is treated as a classical rotator with the 
normalization | S |= 1 and characterized by the angles 
Oi and (pi. Then, semiclassically, the effective hopping 
Hamiltonian is 
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Here Gf creates an electron at site i with spin parallel 
to Si, denotes the nearest-neighbor pairs, which are 
only counted once, and the hopping amplitude acquires 
a Berry's phase and it becomes a complex number given 
by @, 



ti.j = t ( cos — cos -4- 
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This complex hopping appears after rotating the conduc- 
tion electron spins so that the spin quantization axis at 
site i is parallel to Si, and then project onto the spin par- 
allel to Sj. In this paper we will call complex DE (cDE) 
model the system governed by the Hamiltonian ([!]) with 
the hopping given by (g). 

For topologies where electrons close paths do not oc- 
cur, it is possible to choose wave functions phase factors 
such that 



tj i \= t cos ■ 



(3) 



being 9ij the angle between the semiclassical spins Si 
and S,-. 
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In this work we will use the name DE model for the sys- 
tem controlled by equation (|l|) with the hopping (|J) . This 
is the Hamiltonian proposed by Anderson and Hasegawa 
j|| as a generalization of the Hamiltonian describing the 
tunneling between two fixed sites i and j. 

De Gennes [|| and Kubo and Ohata Q| calculated, in 
a mean field theory approximation, the magnetic phase 
diagram of the DE model. Their results show a ferro- 
paramagnetic phase transition at a critical temperature 



MF 



(4) 



This transition is accompanied by a change in the tem- 
perature dependence of the resistivity. In equation (|4|), 
(O)o means the expectation value of the operator O at 
zero temperature. 

However, recently it has been pointed out |8) that the 
critical temperatures predicted by the mean field theory 
of the DE model is much bigger than the observed ex- 
perimentally. Furthermore, the resistivity implied by the 
hopping ([|) is incompatible with many aspects of the ex- 
perimental information. Millis et al |q] proposed that in 
addition to the DE mechanism a strong electron-phonon 
interaction plays a crucial role in these perovskites. The 
interplay between these two effects could reproduce the 
experimental critical temperature and the observed be- 
havior of the resistivity at the magnetic transition tem- 
perature §,0. 

On the other hand, Yunoki et al |ll[ claim that the 
inclusion of the phase in the hopping term can lower 
the critical temperature down to the experimental val- 
ues. Also, for large but finite Hund's coupling, they ob- 
tain that at low concentrations phase separation between 
hole-poor antiferromagnetic and hole-rich ferromagnetic 
regions occurs. The existence of phase separation is also 
obtained by solving the DE model plus an antiferromag- 
netic coupling between next neighbors core spins 

In this paper we are interested in two points; first, in 
comparing the magnetic phase diagram of the two DE 
models, equations (||) and (||), and second, in obtaining 
an approximated expression for the DE energy which al- 
lows us to perform Monte Carlo calculations in bigger 
unit cells and to obtain accurate values of the magnetic 
critical temperature. The main results of our calculations 
are shown in Fig.|l|. Here we show the critical tempera- 
ture versus concentration, for the complex DE model and 
for the real DE model, in the case of a unit cell of size 
4x4x4. We obtain that there is not big difference be- 
tween the critical temperatures of both models. This re- 
sult indicates that the contribution to the partition func- 
tion of electronic configurations in which fermions move 
on closed loops in real space can be neglected. In Fig.[j] 
we also show Monte Carlo results obtained by using an 
approximated expression for the electron kinetic energy. 
We have checked in small systems that this approxima- 
tion is realistic. From our results we obtain that the 
critical temperature of the DE model is in the range of 
the experimental one Pj. 



The paper is organized as follows: in section II we de- 
scribe briefly the Monte Carlo algorithm we use and in 
section III we present the results obtained for the cDE 
and DE Hamiltonians. In section IV we introduce an 
approximation to the DE energy which allows us to per- 
form Monte Carlo simulations in bigger systems. Section 
V is devoted to show the results obtained in the frame of 
this approximation. Finally, we present the conclusions 
in section VI. 



II. DESCRIPTION OF THE MONTE CARLO 
ALGORITHM 

For calculating the magnetic phase diagrams, we have 
performed classical Monte Carlo (MC) simulations on the 
classical core spin angles. The simulations are done in 
N x N x N cubic lattices with periodic boundary condi- 
tions. Although the localized spins are considered clas- 
sical, the kinetic energy of the conduction electrons is 
calculated by diagonalizing the DE Hamiltonian because 
we consider it as a quantum quantity. 

The standard Metropolis algorithm [|l3| was used in the 
MC simulations. The sites to be considered for a change 
in the spin orientation are randomly chosen. Once a site 
is selected for a spin reorientation, the angle associated 
with an attempted change of the spin is chosen at ran- 
dom from within a specified range ]l4| ] . Then the energy 
change, AE, associated with the attempted update, is 
calculated. If the quantity exp{—AE/T) is smaller than 
a random number between and 1, the change is allowed, 
otherwise, it is rejected. Typically, 5000-7000 MC steps 
per spin are used for equilibration and 3000-5000 steps 
for spin are used for calculating averages. 

In the simulations we calculate the average of the in- 
ternal energy, E, and the average of the absolute value 
of the magnetization, M, 
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(5) 
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where ( ) denotes statistical average. Since the MC up- 
dating procedure generates uniform rotations of the spin 
system, a calculation of the MC average of the direction 
of the magnetization is not meaningful. Also we calculate 
the average value of the width of the density of states, 
W, and the value of the chemical potential, /i, measured 
with respect to the bottom of the density of states, 



W=(e 
/' 
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(7) 

(8) 



where e m m and e max are the minimum and the maxi- 
mum energies obtained from diagonalizing the electron 
Hamiltonian, and e OC c is the higher energy of the occu- 
pied states. 



2 



III. MONTE CARLO RESULTS FOR THE 
DOUBLE EXCHANGE MODELS 

In order to obtain the internal energy in each of the 
Monte Carlo steps, it is necessary to diagonalize the elec- 
tronic Hamiltonian. The diagonalization is very expen- 
sive in terms of CPU time, and therefore we can only 
study small size systems. Here we present the results 
obtained for a cubic unit cell of size N = 4 and for the 
electron concentrations x = 0.1, 0.25 and 0.3. 

In figure ||, it is presented, for the DE and the cDE 
Hamiltonians, the magnetization, M, as a function of 
the temperature, T. Due to the finite size of the unit 
cells used in the simulations, M is different from zero at 
any temperature, and we define the critical temperature, 
T c , as the point where the second derivative of M with 
respect T changes sign. This way of obtaining the crit- 
ical temperatures, implies uncertainties of around 10% 
of the value of T c . In any case we obtain that the crit- 
ical temperatures of the DE and the cDE Hamiltonians 
are practically the same. In figure |^ we plot the internal 
energy E as a function of T. Note that the difference 
between the internal energies obtained by using the cDE 
and using the DE Hamiltonians is rather small. That is 
the reason why both models give very similar magnetic 
critical temperatures. 

These results imply that in order to calculate total 
kinetic energies the close paths are almost irrelevant. 

In figure ^ we plot, as a function of the temperature, 
the width of the density of states W and the chemical 
potential fx, for the DE case. Note that at temperatures 
close to T c , the bandwidth and the chemical potential 
have a value bigger than the obtained in the fully disor- 
dered case, T — > oo. This is because in the DE models 
the ferromagnetic to paramagnetic transition is a second 
order phase transition, and therefore the internal energy 
is continuous at the critical temperature. Moreover, the 
bandwidth and the chemical potential present a contin- 
uous behavior near T c . In the DE models we do not 
observe any change in the electronic states at the ferro- 
paramagnetic critical temperature. 



IV. APPROXIMATION FOR THE DOUBLE 
EXCHANGE ENERGY 

In the MC calculations of the DE models, the size of 
the matrix to diagonalize impose a restriction on the di- 
mension of the unit cell used in the simulation. In order 
to be able to perform simulations in bigger systems, we 
have developed a second order perturbation theory for 
obtaining an expression for the electron kinetic energy in 
a background of randomly oriented core spins. 

We start writing the Hamiltonian (111) in the form, 



H = H + V 



(9) 



H = (C+C 3 + h.c.) 

V^-^2{St i>j C+C j + h.c. 

(id) 



(10) 
(11) 



being t the average of the absolute value of the hopping 
amplitude, 



t 



3N, 



and 



(12) 



(13) 



Here Nq is the total number of Mn ions in the system. 
Note that since the Hamiltonian is hermitic, t is a real 
quantity. 

Given a disordered system, characterized by a set of 
{tij}, we want to obtain the expectation value of H. 
The Hamiltonian Hq can be diagonalized by using the 
Bloch states, 



k) = 




(14) 



where | i) represents the atomic orbital at site i, Ri are 
the lattice vectors and k is a wave vector in the first 
Brillouin zone. The energy of the state | k) is 



e(k) = — 2t (cos k x a + cos k y a + cos k z a) 



(15) 



being a the lattice parameter. In function of these eigen- 
values, the expectation value of Hq is, 



(16) 



where the sum is over the occupied states. From the value 
of Eq we can obtain the expectation value of C^Cj, 



(C+Q) = - 



Eq 
6iN ' 



(17) 



which is independent of the value of i. 

To obtain the expectation value of (V), we first notice 
that since the set {SUj} are randomly distributed, (k | 
V | k) = 0, and it is necessary second order perturbation 
theory in order to get a correction to the expectation 
value of H , 



E 2 ~-J2{k\VG (e k )V\k) , 



where 



with 



G (M = 
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is the Green function of the perfect crystal. Because the 
values {Stij} are not correlated, it is easy to obtain the 
expression, 

^2 = -a 2 ^(i-t 2J ) 2 , (20) 

(id) 

with 

1 , e(k + k') \ 1 

(21) 

Adding Eq and E%, we obtain the following expression 
for the energy of the system, 

S^-2<C+Q)o^cos^-a 2 ^(?-^) 2 . (22) 

(id) (id) 

In figure || we plot (C^Cj)o and a 2 as a function of the 
electron density x. The first term in the above expression 
is a ferromagnetic coupling for the core spins and the 
second term is an antiferromagnetic coupling. Note that 
(CfCj)o is only linear on x at small values of the electron 
concentration. 



V. MONTE CARLO RESULTS WITH THE 
APPROXIMATED DOUBLE EXCHANGE 
ENERGY 

In figure || we plot, for x — 0.1, x — 0.25 and x — 0.3, 
the absolute value of the magnetization, M, as a function 
of the temperature, T, for the case where the internal en- 
ergy is obtained by using the approximated DE energy, 
equation (|2|) , for comparing with the case in which it is 
obtained by diagonalizing the DE Hamiltonian. The dif- 
ference between both curves is very small, and the corre- 
sponding critical temperatures are practically the same, 
see figure [lj. Therefore we conclude that equation is 
a good approximation for the energy of the DE model. 

Using equation (^2|) , for the DE energy, we can perform 
MC simulations in much bigger systems, and obtain more 
precise values of the critical temperatures. In figure ^ we 
plot the magnetization as a function of the temperature 
for a sample of size N = 20. For this size the critical 
temperature is obtained with a big precision. We have 
checked that for N = 20 the values of the critical tem- 
peratures are accurate in two digits and that it is not 
worthwhile to increase more the size of the system. We 
have performed MC simulations for different values of the 
electron concentration and, in figure [j], we have plotted 
the magnetic phase diagram. 

Kubo and Ohata have developed an expression for 
the DE energy which coincides with Eq . In order to know 
the importance of E^ in the calculation of the value of the 
critical temperature, we have also computed the value of 



M as a function of T by using only Eq for the internal 
energy, see figure |[ Comparing these results with the 
ones obtained with the approximated DE energy, E0+E2, 
we estimate that the inclusion of E2 in the internal energy 
lowers the critical temperature more than a 20%. 

In Fig.^, we plot for different values of the concentra- 
tion x, the statistical average of the absolute value of 
the hopping amplitude, as a function of the tem- 

perature. This quantity is proportional to the electron 
bandwidth. We find that is a continuous func- 

tion of T, and there are not signals of any change in the 
electronic structure near T c . 

The bandwidth near T c is around 1.15 times bi gger 
than the bandwidth in the T — > 00 limit. Li et al |lq | 
have obtained that in the T — » 00 limit of the cDE model, 
only ~ 0.5% of the electron states are localized. In order 
to know the difference between the density of states at 
T r*> T c and at T — > 00, we have calculated the statistical 
average of the electron density of states. We find that the 
number of electron states with energies between —W(T C ) 
and —W(T — » 00) is less than 1% of the total number of 
states. Therefore, even if all these states were localized, 
the percentage of localized states near T c is not larger 
than 2%. We conclude that in the DE models it is not 
possible to relate the metal-insulator transition with the 
ferro-paramagnetic one. 

The critical temperatures we obtain are around 1.5 
times smaller than the obtained in mean field theory 
by Kubo and Ohata [j]]. However, it is around 8 times 
smaller than the obtained by Millis et al. || using also 
a mean field approximation. We think that the differ- 
ence occurs because they use a lineal dependence, see 
Fig. ||, of the kinetic electron energy which overestimate 
the value of (Cj~Cj)o, and therefore the value of T c . 

From Kubo and Ohata model ([|), it is clear that T c 
scales with t(Cj~Cj)o and the magnetization versus T/T c 
curve is independent of the electron concentration x. But 
there is no reason to think that this must be the case 
when the second term, E2, of our approach is included 
in the calculations. However, we have found that the M 
versus T/T c curves almost coincide (see Fig.^J). This fact 
is not only found in the frame of our approximation to the 
DE energy but also when the DE and cDE Hamiltonians 
are diagonalized. 

The presented results are in qualitatively agreement 
with those obtained by Dagotto et al (T^|. However our 
estimations of T c are ~ 1.5 times bigger. The discrepancy 
is due to the difference in the criterion used to obtain T c . 
In reference Jl(| the T c is defined as the temperature 
where the spin-spin correlation function in real space be- 
comes zero at the maximum distance available in the unit 
cell. They used a unit cell of size 6x6x6. Our criterium 
is based in the change of sign of the second derivative 
of M with respect to T. We have checked that both 
criteria give the same T c in big unit cells, but the cri- 
terium based on the correlation function underestimate 
T c in small unit cell calculations. This is clear in Fig.|§| 
where we plot the spin-spin correlation function versus 
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distance for several temperatures and for two unit cells 
with TV = 6 and N = 10. The electron concentration is 
x = 0.14 to compare with the results presented in refer- 
ence Jl6| . In the case of N = 6 it seems that the critical 
temperature is smaller than T = 1/15, because the cor- 
relation function is zero at a separation of 8a. However, 
when N — 10 it is clear that the ferromagnetic correla- 
tion at T = 1/15 is not zero at 8a, and in fact it seems 
that at this temperature the system is still ferromagnetic. 
On the other hand, our criterium overestimate the T c for 
small unit cells and the possibility of increasing the size 
leads to more realistic values of the critical temperature. 

To compare with the experimental critical tempera- 
tures, it is necessary an estimation of t. The experimental 
bandwidth is around l-4eV 0, and therefore the value 
of t is around 0.08-0.3eV. With this result our estimate 
of the critical temperature is between 150K and 500K at 
x ~ 0.2. This value of T c is in the range of the observed 
experimentally. 



VI. CONCLUSIONS 

Using Monte Carlo techniques, we have obtained the 
magnetic phase diagram of the double exchange Hamil- 
tonian. Comparing the results obtained from the double 
exchange Hamiltonian with a complex hopping (cDE) 
and with its absolute value (DE) we have found that 
the Berry phase of the hopping has a negligible effect in 
the magnetic critical temperature. This implies that it 
is possible to choose wave functions phase factors that 
counteract the phase of the hopping such that it is not 
relevant in the MC calculations. To avoid the limitations 
on the size of the systems studied we have developed a 
second order perturbative approach to calculate the elec- 
tron kinetic energy without diagonalizing the Hamilto- 
nian. Within this approach we have calculated the criti- 
cal temperature for a bigger size of the system (N = 20) 
and, consequently, more accurately. The values of T c ob- 
tained are in the range of the experimental ones and in 
qualitative agreement with those in ref . |l6| . Neither the 
average of the electron density of states nor the kinetic 
energy reveal any strange behavior near the magnetic 
critical temperature. We conclude that, in order to un- 
derstand the metal-insulator transition which occurs in 
the manganites, another effects, not included in the dou- 
ble exchange Hamiltonian, should be taken into account. 
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FIG. 1. Critical temperature of the ferro-paramagnetic 
transition versus concentration of conduction electrons for dif- 
ferent approaches. Squares correspond to the Hamiltonian (|l]) 
with the hopping given by (^); circles to hopping Q. The 
curves are related to the calculations made with an approxi- 
mated expression obtained in order to avoid diagonalization 
in the MC simulations (ADEE: approximation to the double 
exchange energy). This fact allows us to increase the size of 
the system. The line is only a guide to the eye. The critical 
temperatures obtained are in the range of the experimental 



FIG. 2. Absolute value of the magnetization versus the 
temperature for concentrations x = 0.1, 0.25 and 0.3. The 
size of the unit cell is TV = 4. In order to compare we have 
plotted also the curves obtained with the cDE and DE Hamil- 
tonians as well as with the perturbative approach developed 
here (ADEE). 



FIG. 3. Double exchange energy versus temperature. The 
lack of difference in these curves implies the result in Fig. ^. 
The corresponding critical temperatures are pointed with an 
arrow. 
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FIG. 4. Bandwidth (a) and chemical potential (b) versus 
temperature. In the disordered case, T — * oo, W should be 
equal to 8. Here, in the paramagnetic phase, W > 8 because 
the transition is second order (see text for details). Note that 
the curves are continuous near T c which is pointed with an 
arrow in each case. 



FIG. 5. Coefficients used in the perturbative approach to 
the DE Hamiltonian versus concentration (see eq. (2^)). 



FIG. 6. Magnetization as a function of T/T c as obtained 
using the th order approximation (dashed line), Eq in eq. 
(pj[), and the 2 nd order approximation (continuous line) to 
the double exchange energy (ADEE), eq. (p^). The critical 
temperature taken for reference is the one that corresponds 
to the ADEE model. In the case of the th order approxima- 
tion, the M versus T/T c curve is independent of the electron 
concentration. In the ADEE model, there is no reason for this 
independence; however, we have found that M(T/T C ) is prac- 
tically independent of x. To show this, the M(T/T C ) curve 
has been plotted for various concentrations (x = 0.1, 0.2, 0.3, 
0.4 and 0.5). 



FIG. 7. Average of the absolute value of the hopping am- 
plitude, versus temperature for N = 20 and concen- 
trations x — 0.1, 0.25 and 0.3. This average is proportional 
to the bandwidth. We recover the continuous behavior that 
we found on Fig. 4. 



FIG. 8. The spin-spin correlation versus distance is plot- 
ted for different temperatures. In (a) N = 6 and in (b) 
N — 10. For these calculations we have used the perturbative 
approach. If the criterium to obtain T c is that the spin-spin 
correlation becomes zero, a smaller size of the system leads 
to an underestimation of the critical temperature. 
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Fig. 1 M.J. Calderon and L. Brey 




Fig. 3 M.J. Calderon and L. Brey 
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Fig. 6 M.C. Calderon and L. Brey 
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